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Abstract. Sensitivities of the oceanic biological pump within 
the GISS (Goddard Institute for Space Studies ) climate mod- 
eling system are explored here. Results are presented from 
twin control simulations of the air-sea CO 2 gas exchange 
using two different ocean models coupled to the same atmo- 
sphere. The two ocean models (Russell ocean model and Hy- 
brid Coordinate Ocean Model, HYCOM) use different ver- 
tical coordinate systems, and therefore different representa- 
tions of column physics. Both variants of the GISS climate 
model are coupled to the same ocean biogeochemistry mod- 
ule (the NASA Ocean Biogeochemistry Model, NOBM), 
which computes prognostic distributions for biotic and abi- 
otic fields that influence the air-sea flux of CO 2 and the deep 
ocean carbon transport and storage. In particular, the model 
differences due to remineralization rate changes are com- 
pared to differences attributed to physical processes modeled 
differently in the two ocean models such as ventilation, mix- 
ing, eddy stirring and vertical advection. GISSEH(GISSER) 
is found to underestimate mixed layer depth compared to 
observations by about 55% (10%) in the Southern Ocean 
and overestimate it by about 17% (underestimate by 2%) 
in the northern high latitudes. Everywhere else in the global 
ocean, the two models underestimate the surface mixing by 
about 12-34 %, which prevents deep nutrients from reaching 
the surface and promoting primary production there. Conse- 
quently, carbon export is reduced because of reduced produc- 
tion at the surface. Furthermore, carbon export is particularly 
sensitive to remineralization rate changes in the frontal re- 
gions of the subtropical gyres and at the Equator and this 
sensitivity in the model is much higher than the sensitivity to 
physical processes such as vertical mixing, vertical advection 


and mesoscale eddy transport. At depth, GISSER, which has 
a significant warm bias, remineralizes nutrients and carbon 
faster thereby producing more nutrients and carbon at depth, 
which eventually resurfaces with the global thermohaline cir- 
culation especially in the Southern Ocean. Because of the 
reduced primary production and carbon export in GISSEH 
compared to GISSER, the biological pump efficiency, i.e., 
the ratio of primary production and carbon export at 75 m, is 
half in the GISSEH of that in GISSER, The Southern Ocean 
emerges as a key region where the CO 2 flux is as sensitive to 
biological parameterizations as it is to physical parameteri- 
zations. The fidelity of ocean mixing in the Southern Ocean 
compared to observations is shown to be a good indicator of 
the magnitude of the biological pump efficiency regardless 
of physical model choice. 


1 Introduction 

The world’s oceans help regulate the amount of heat, CO 2 
and other greenhouse gases in the atmosphere that are re- 
sponsible for the radiative imbalance of the planet. Oceans 
therefore play a crucial role in the global climate. There 
are two processes that control the air-sea exchanges of 
CO 2 : the solubility pump, and the biological (and carbon- 
ate) pump (Volk and Hoffert, 1985; Gruber and Sarmiento, 
2002; Falkowski et ah, 2003). About one-third of the an- 
thropogenic long-term averaged CO 2 emissions released into 
the atmosphere is absorbed by the ocean (Khatiwala et al. 
(2009) and references therein) whereas shorter-term aver- 
ages amount to one-fourth (Global Carbon Project, 2012; 
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LeQuere et al., 2012) through the solubility pump. The net 
uptake flux of CO 2 by the global oceans is estimated to be 
about 1.6±0.9PgCyr“* (Takahashi et ah, 2009) with the 
Southern Hemisphere oceans being the largest CO 2 sink tak- 
ing up about 1 . 1 Pg Cyr“ * , while the northern subtropical and 
midlatitude oceans take up about 0.7PgCyr“*. The equato- 
rial oceans emit 0.7PgCyr“* to the atmosphere. Uncertain- 
ties in the solubility pump stem mainly from uncertainty in 
the pC 02 distributions at the surface and the parameteriza- 
tion of the gas exchange at the air-sea surface. These range 
from 3 to 30% with signiflcant regional differences (Sig- 
norini and Mcclain, 2009; Feely et al., 2004). 

There is greater uncertainty in our knowledge of the ocean 
biological pump (OBP), which describes the amount of car- 
bon flxed by biology at the surface and then redistributed 
within the ocean (Denman, 2003; Buesseler et al., 2008; 
Henson et al., 2011), i.e., those processes, such as detri- 
tal sinking and remineralization, that remove surface carbon 
from the euphotic zone of the ocean and transport it to the 
deeper layers. It is estimated that this deep ocean carbon ex- 
port is about 6PgCyr“^ Martin et al. (1987) to 10PgCyr“' 
(Falkowski et al., 1998) usually depending on the type of ob- 
servations while model results range from about 8PgCyr“^ 
to llPgCyr“^ (Taucher and Oschlies, 2011; Dunne et al., 
2007) depending on the model parameterizations. 

To be able to reduce the uncertainty range and the biases 
in the model results with respect to the OBP, we need to 
understand the type of sensitivities different parameteriza- 
tions of the OBP produce, and how these sensitivities com- 
pare to other model biases, e.g., due to model physics and 
how much of these sensitivities are model dependent. Sev- 
eral recent studies (Oschlies, 2001; Carr et al., 2006; Kwon 
et al., 2009; Taucher and Oschlies, 2011; Kriest et al., 2010) 
showed that the biologial pump has important sensitivities 
to model parameters such as temperature, remineralization 
depth, nutrient uptake/growth rate and sinking of particu- 
late matter. These key process parameterizations directly af- 
fect the biological pump, and through changes in the surface 
pC 02 distributions, they have indirect effects on the solubil- 
ity pump. 

Studies based on simple biogeochemical models with few 
components (Bacastow and Maier-Reimer, 1991; Kwon and 
Primeau, 2008) demonstrated that when the remineralization 
length scale is increased, or equivalently the remineralization 
rate is decreased, near the surface there are increases in the 
meso- and bathypelagic nutrient distributions and the nutri- 
ent transport, which traces the large-scale deep circulation 
pathways. 

More complex biogeochemistry model studies showed 
that even small changes in remineralization depth lead to re- 
distribution of recycled carbon from the meso- to the bathy- 
pelagic waters and lead to reduction of the atmospheric CO 2 
(Kwon et al., 2009). 

Preformed nutrients are remineralized in the epipelagic 
zone within the young North Atlantic Deep Water (NADW) 
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and subsequently stored in the oldest waters near the ocean 
bottom. Therefore, the NADW is a major pathway of atmo- 
spheric CO 2 being sequestered into the deep ocean. 

Kriest et al. (2012) showed that global biases in observed 
nutrient distributions depend mainly on parameter values of 
speciflc processes related to recycling and sinking of nu- 
trients, rather than on other processes in the biogeochem- 
istry model. In their experiments, simple parameterizations 
of the remineralization and sinking led to large global biases 
from observations due to deflcits in nutrient concentrations 
in the intermediate depths of the ocean due to excessively 
large export of organic matter to greater depths. Moreover, 
they found that even though global biases are low, this may 
not represent the overall skill of the parameterization, as re- 
gional biases often cancel each other. North Pacific deep wa- 
ters are very old and therefore all biases, related to biology 
and physics model deficiencies, slowly accumulate there over 
time (Khatiwala et al., 2012). Changes in these parameteri- 
zations account for most of the model biases in the North Pa- 
cific. Therefore, the simpler models, models with fewer com- 
ponents, have better skill in the representation of the surface 
properties in this region, at depth nonlocal errors (i.e., model 
errors from remote regions) are larger. Kriest et al. (2012) 
also And that adjusting these parameters may enhance phy- 
toplankton production particularly in the high nutrient low 
chlorophyll (HNLC) regions (e.g., the Southern Ocean, the 
equatorial regions or the northern North Pacific regions) but 
this may compensate for the lack of other nutrient limitations 
in the model (e.g., iron). 

Several studies (Kriest et al., 2010; Kriest and Oschlies, 
2011) have drawn attention to the role of physical process 
parameterizations in the models such as vertical mixing, or 
even numerical diffusion by the model advection scheme, 
(e.g., the upstream advection operator) in transporting some 
of the particulate organic matter to depth, which may either 
compensate or enhance biases due to biological parameter- 
izations. Similarly, low model vertical resolution leads to 
increased particle flux via increased numerical mixing, and 
therefore leads to more recycling at depth, compared to a 
model with a finer vertical resolution. 

Finally, improving the model fidelity to observed nutrient 
distributions by adjusting the parameterizations of the verti- 
cal particle flux may degrade simulations of other fields such 
as primary production (PP) (Carr et al., 2006; Kriest et al., 
2010). In these studies, it is shown that a model that predicts 
PP skillfully may not perform as well in simulating nutrient 
concentrations. 

In this paper we will focus on one biogeochemistry model 
and vary the remineralization rate within two different ocean 
models coupled to the same atmosphere. Since we are able 
to use two ocean models that represent different classes of 
ocean modeling (z coordinate vs isopycnal) we are able to 
assess the relative importance of the biological pump sensi- 
tivities to biological and physical parameterizations. 
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The paper is structured as follows: Sect. 2 includes infor- 
mation on the GISS climate model configurations with the 
two ocean models and the biogeochemistry module. Results 
from two sets of twin runs where we vary the remineral- 
ization rate are presented in Sect. 3. We distinguish effects 
of the remineralization rate change on the biological pump 
(Sect. 3.1), i.e., detritus, nitrates, total chlorophyll, primary 
production fields and deep carbon export, and effects on the 
gas exchange pump (Sect. 3.2), i.e., the dissolved inorganic 
carbon (DIG), pC02 distributions and the air-sea flux of 
CO 2 . Our results are summarized and are placed in the con- 
text of similar studies in the Discussion section (Sect. 4) and 
lastly some conclusions are offered in Sect. 5. 


2 Model description 

The carbon cycle runs are described extensively in Romanou 
et al. (2013) and have been submitted to the Coupled Model 
Intercomparison Project-Phase 5 (CMIP5) as GISS-E2-R- 
CC and GISS-E2-H-CC, which are the carbon cycle (CC) 
versions of model E2 coupled to the Russell (R) or HYCOM 
(H) ocean models. All carbon cycle runs started from equi- 
librium conditions of coupled atmosphere-ocean-ice 3000 yr 
runs, and were subsequently integrated again with fully prog- 
nostic carbon cycle to equilibrium for an additional 300 yr. 
Here we present only preindustrial control CO 2 simulations 
where the atmospheric pC 02 responds to the ocean, but the 
global average is fixed to 285.2 ppmv (parts per million by 
volume). A very brief description of the two ocean and the 
biogeochemical models is given here. 

The Russell ocean model is a non-Boussinesq, mass- 
conserving ocean model with up to 32 vertical levels and 
1° X 1.25° horizontal resolution. The vertical coordinate is 
mass per unit area because mass is formally conserved in 
the model and therefore introduction of freshwater fluxes 
can be easily handled. Normalizing to grid box area ensures 
that computations and comparisons of fluxes to/from adja- 
cent boxes are properly handled. The model has a free sur- 
face and natural surface boundary fluxes of freshwater and 
heat whose prognostic variable is potential enthalpy. It uses 
the Gent and McWilliams (1990) scheme for isopycnal eddy 
fluxes and isopycnal thickness diffusion. Vertical mixing is 
done according to the KPP (K-Proflle Parameterization) ver- 
tical mixing scheme (Large et al., 1994) with nonlocal mix- 
ing effects included. Vertical advection is computed using a 
centered differencing scheme. 

The Hybrid Coordinate Ocean Model (HYCOM, Sun and 
Bleck (2006) variant C) uses pressure coordinates near the 
surface and isopycnals in the ocean interior. Horizontal reso- 
lution is 1° X l°cos(0), where 4> is the latitude, vertical reso- 
lution of 26 coordinate layers. The vertical mixing scheme 
is KPP (Large et al., 1994) near the surface and follow- 
ing McDougall and Dewar (1998) in the isopycnic interior. 
Deep convection at high latitudes is parameterized through a 
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convective adjustment scheme. Brine rejected during freez- 
ing is uniformly distributed in the uppermost 200 m when 
water depth is greater than 1000 m or the total depth oth- 
erwise. The coefficient of lateral T /S diffusion is increased 
near the Equator to account for the tropical instability waves, 
which are very diffusive and which are not resolved by the 
coarse grid. The Gent-Mc Williams parameterization is im- 
plemented via interface smoothing, and its effect therefore is 
limited to the isopycnic coordinate domain (which excludes 
the near-surface waters). In the present implementation, only 
mass fluxes at even time steps are used in building up the 
mass flux time integrals. This is to avoid leapfrog-related in- 
consistencies between layer thickness tendencies and hori- 
zontal flux fields, which would generate spurious diapycnal 
mass fluxes. Their sum is then used at the end of a longer time 
step (1/2 day) to rigorously reconstruct all terms of the time- 
integrated continuity equation. (The other terms are the net 
change in layer thickness over the long time step and the ver- 
tical mass fluxes that are derived diagnostically from the oth- 
ers.) This allows us to solve the flux form of the tracer advec- 
tion equation over a large time step while rigorously conserv- 
ing everything in an environment of perpetually shifting layer 
configurations. The single remaining inconsistency is the one 
caused by time smoothing of the thickness held. It spawns 
small diapycnal fluxes (and hence interlayer tracer leakage) 
even in the absence of physically based vertical mixing pro- 
cesses. The numerical method used for 3-D tracer advection 
is FCT (flux corrected transport), which maintains positive- 
definiteness and monotonicity in the tracer fields. The hori- 
zontal fluxes, before being clipped, are 2nd order centered, 
and the vertical fluxes (again before clipping) are computed 
using the piecewise parabolic method (PPM). Solving trans- 
port equations in flux form does not rigorously conserve total 
tracer when layer thickness is allowed to go to zero. Conser- 
vation errors are currently eliminated by applying a global 
correction factor to each tracer held following each transport 
step. With this corrective device in place, tracer is rigorously 
conserved in the model. 

Most relevant to the present study is the fact that the differ- 
ent coordinate system in the two ocean models produces dif- 
ferent distributions of seawater properties below the mixed 
layer depth where numerical and spurious diapycnal mix- 
ing can become important (Bleck et al., 1992; Griffies et al., 
2000). Other relevant differences between the two ocean 
models are the extent of implementation of the vertical mix- 
ing scheme (KPP); in GISSER KPP extends to the bottom 
of the water column but nonlocal terms are not included, 
whereas in GISSEH it only reaches the bottom of the mixed 
layer. Eddy mixing along isopycnals (Gent and McWilliams, 
1990) is only explicitly modeled in GISSER, whereas in GIS- 
SEH it is approximated by isopycnal thickness diffusion for 
the regions where these isopycnals do not outcrop. More- 
over, there is no tracer horizontal advection in GISSEH. Al- 
though the models are different in other ways, the signifi- 
cantly decreased vertical diffusion/mixing below the mixed 
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layer in the GISSEH is perhaps the main difference between 
them. 

The interactive ocean carbon cycle model is the NASA 
Ocean Biogeochemistry Model (NOBM, (Gregg and Casey, 
2007; Romanou et ah, 2013), which is interactively coupled 
to the ocean model from which it uses temperature, salin- 
ity, mixed layer depth, the ocean circulation fields, horizon- 
tal advection and vertical mixing schemes and the atmo- 
spheric model from which it uses the surface shortwave ra- 
diation (direct and diffuse) and surface wind speed. NOBM 
includes four phytoplankton groups (diatoms, chlorophytes, 
cyanobacteria, and coccolithophores), four nutrients (ni- 
trate, silicate, ammonium and iron), three detrital pools (ni- 
trates/carbon, silicate and iron) and a single herbivore com- 
ponent. The carbon submodel includes dissolved organic and 
inorganic carbon (DOC and DIC, respectively) and surface 
alkalinity, which is a constant function of salinity at the sur- 
face. The underwater irradiance needed to drive phytoplank- 
ton growth is calculated explicitly for the entire spectrum of 
light with the use of a radiative transfer model that takes into 
account water absorption and scattering based on the opti- 
cal properties of the phytoplankton groups that were derived 
from laboratory studies (Gregg, 2002). Detrital settling and 
phytoplankton sinking are modeled in HYCOM using a one- 
dimensional piecewise parabolic transport method whereas 
in the Russell ocean model using a piecewise linear method. 

In Gregg and Casey (2007) phytoplankton sinking was 
modeled as an additional advection in the vertical. Sink- 
ing rates were specified at 31 °C and derived from Stokes’ 
law using representative phytoplankton sizes for each group 
from a very thorough survey of the literature (Gregg and 
Casey (2007) and references therein). The sinking rates Wg 
of each phytoplankton group were adjusted to account for 
viscous effects, according to Stokes’ law as a function of 
temperature T : 

Ws(r) = w;s(31)[0.451-h0.0178 J]. (1) 

Coccolithophores were assumed to be able to sink at speeds 
that were a function of the growth rate from 0.3 to 1 .4 m d“ ^ 
based on observations by Fritz and Balch (1996). A linear 
relationship was assumed: 

w^seocci = 0.752p^,o^,^,i (high) -F 0.225, (2) 

where Wg is the sinking rate of coccolithophores (m d“*) and 
p(high) is the maximum growth rate actually achieved for the 
previous day. 

For three of the nutrient groups, i.e., nitrate, silica, and 
dissolved iron, corresponding detrital pools are defined in or- 
der to simulate the remineralization and return to dissolved 
form. The equations for these detrital pools are given in the 
Appendix of Gregg and Casey (2007), Eqs. (A7-A9). Only 
carbon detritus is explicitly modeled in NOBM since the 
C : N ratio is constant and conversion to the nitrate pool is 
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straight forward. Detritus pools are changed due to advec- 
tion, diffusion and sinking. Similar to phytoplankton sink- 
ing, detrital sinking depends on viscosity and is a func- 
tion of temperature. Remineralization is also temperature- 
dependent and varies with the phytoplankton growth rate for 
each phytoplankton group. The remineralization rate for ni- 
trogen/carbon is 0.1 d“*, for silica is 0.002 d“' and for iron 
is0.7d“'. 

3 Results 

Results are presented below for two sets of model runs for 
each climate model configuration. GISSER-l(h) stands for 
Russell ocean model with low(high) remineralization rate 
and GISSEH-l(h) stands for HYCOM and low(high) rem- 
ineralization rate. 

Remineralization rates in the ocean are highly uncertain 
and vary significantly regionally as the VERTIGO experi- 
ment has demonstrated (Buesseler et al., 2008). Tian et al. 
(2000) showed that regionally (northern Atlantic subtropi- 
cal gyre) the remineralization rate in the upper 400 m can be 
as high as 0.8d“*. Christodoulaki et al. (2013) suggested a 
value of 0.45 d“*. At the same time, little is known about 
how and how much such large regional differences affect es- 
timates of the biological pump and the air-sea gas exchange. 
We chose here to present results for two cases: a “low” rem- 
ineralization rate of 0.1 d“*, which is the value suggested 
by Denman (2003), and a “high” rate of 0.5 d“^ Analysis 
is extended to cover more cases in order to investigate the ro- 
bustness of our conclusions, and the parameter range extends 
from 0.01 to 0.5 d“^ 

Both NASA-GISS ocean-atmosphere climate models are 
coupled to NOBM and run out to equilibrium under 
constant preindustrial atmospheric CO 2 concentrations of 
285.23 ppmv. Equilibrium is achieved when the gas ex- 
change between the ocean and the atmosphere stabilizes to 
a value about 0.9 Pg Cyr“^ . The differences in the biological 
and gas exchange pumps between the two models are then 
evaluated. Detailed model assessment is given in Romanou 
et al. (2013). Here we will focus on the differences between 
the two physical models, particularly in the regions/cases 
where there are large biological pump differences as well. 

3.1 Sensitivity of the biological pump 

3.1.1 Large-scale mixing 

To cast the biological pump sensitivity in light of the physi- 
cal model differences, we review some of the main findings 
of Romanou et al. (2013) where a more detailed model eval- 
uation was presented. In that paper, the model ocean carbon 
cycle is shown to be very sensitive to biases at the surface 
but also the interior of the ocean. Surface thermal, radiative 
and momentum biases give rise to mixed layer depth biases 
whereas deeper ocean vertical mixing differences explain the 
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Nmidlat Nsubtro Ntropic Stropic Ssubtro Smidlat 


Fig. 1. Annual mean mixed layer depth comparison in the mod- 
els and the observations in different regions (zonal averages). The 
northern midlatitude region (Nmidlat) extends from 40° N to the 
pole, the northern subtropical region (Nsubtro) extends from 15° to 
40° N, the northern tropical region (Ntropic) extends from the Equa- 
tor to 15° N, the southern tropical region (Stropic) extends from the 
Equator to 15° S, the southern subtropical region (Ssubtro) extends 
from 15° to 40° S and the southern midlatitude region (Smidlat) 
extends from 40° S to the Antarctic coast. The observations come 
from the deBoyer atlas (deBoyer Montegut et al., 2004). The units 
of the vertical axis are in meters. 

the lower nutrient load that upwells to the surfaee to sustain 
the primary produetion in the GISSEH model. 

Figure 1 shows a regional comparison and evaluation of 
the model mixed layer depth (MLD). Both models under- 
estimate MLD in the Southern Flemisphere, particularly in 
the Southern Flemisphere ’s midlatitudes, where the tempera- 
ture bias in GISSER is 2 °C and in GISSEH is 3 °C (Fig. 2). 
In Romanou et al. (2013) it was shown that wind speeds 
over the Southern Ocean are underestimated in both models 
(Fig. 2) which, compounded by overestimated heat losses in 
the winter and cloud coverage during summer, lead to shal- 
lower mixed layer depths. 

Below 100 m and down to 3000 m (Fig. 2), GISSER ex- 
hibits significant warm bias whereas GISSEH is in better 
agreement with the observations. The GISSER deep warm 
bias is due to significant eddy transport in the Southern 
Ocean. 

3.1.2 Detritus 

In NOBM, detritus, i.e., the nonliving part of the particulate 
organic carbon (POC) pool, is changed by advection, diffu- 
sion, detrital sinking, remineralization, biological production 
and bacterial degradation according to 

dDt^V(KVD)-V-VD-wdDz-RaD (3) 

-|-4’[/c P; + H] -|- 4>(1 — e)r]2H^ — XD, 
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where, D is the detritus concentration, K is the vertical dif- 
fusivity, V is the horizontal current velocity, w is the vertical 
advection (sinking) speed, R is a temperature dependence 
term, a is the remineralization rate at 20 °C, <1> is the car- 
bon : chlorophyll ratio, k is the senescence rate, P, is the vec- 
tor of phytoplankton species, where i is for diatoms, chloro- 
phytes, cyanobacteria, and coccolithophores, rji and r ]2 are 
heterotrophic loss rates, H is the herbivore concentration, e 
is an excretion rate for nitrate/carbon and X is the detrital 
breakdown rate at 20 °C. The values of the parameters are 
given in Gregg and Casey (2007). 

As we increase the remineralization rate, the sink term, 
RaD, on the right-hand side of Eq. (3) increases, and that 
leads to decreased detritus concentrations at all levels but 
more surface intensified (Fig. 3). Most of the reduction oc- 
curs in the Southern Ocean north of 60° S, within the wind- 
driven divergence zone and up to the subtropical convergence 
zone (STCZ) along 40° S approximately. The distributions 
of detritus in the Southern Ocean upwelling and divergence 
zone reflect the effect of changes in remineralization rate and 
the combined effect of the thermohaline circulation pattern: 
as we increase the remineralization rate in the North Atlantic 
convection regions, less detrital material reaches deeper wa- 
ters, therefore less is remineralized back to inorganic forms 
along the deep water pathways and eventually the NADW 
carries less detritus as it upwells in the Southern Ocean di- 
vergence zone. 

In the Southern Ocean, the GISSEH model shows the 
largest detritus reduction in the region that lies between 
the Antarctic Polar Front (APF) at about 50° S and the 
STCZ (?«40° S, Fig. 3). Similar reductions are shown in the 
GISSER model but over a wider region (70°^0° S). The re- 
duction region in GISSEH is narrower than in GISSER, be- 
cause the former model has significantly shallower mixed 
layer depths in the region 70°-50° S since the temperature 
there is much higher than in GISSER (see Fig. 3 in Romanou 
et al. (2013)). Higher temperatures in GISSEH lead to higher 
remineralization rates and therefore lower detritus distribu- 
tions than in GISSER (see Fig. 17 in Romanou et al. (2013)). 

A close inspection of Eq. (3) shows that in NOBM the 
remineralization sink term scales as the local temperature 
profile, through the Bissinger et al. (2008) parameterization, 
which is a refinement to the Eppley curve (Eppley, 1972) for 
most ocean areas. 

If we define “effective remineralization” as (Zeff = Ra : 

(4) 

where a = 0.1 d“* (the low remineralization case), then at 
the surface (Fig. 4 top panels) the remineralization term is 
lower in the tropical upwelling regions and in the polar re- 
gions where the surface waters are colder than in the rest of 
the ocean. The effective remineralization distributions follow 
the SST (sea surface temperature) climatologies, shown in 
Fig. 2 of Romanou et al. (2013). Fronts in temperature corre- 
spond to fronts in the remineralization rate. The two models 
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Fig. 2. Potential temperature profiles from the two models and observations in the same regions as in Fig. 1. Vertical axis is stretched near 
the surface and the units are in meters. 



Fig. 3. Surface distributions of SDetritus in (a) the GISSER (shorthand: R) and (b) the GISSEH (shorthand: H) models at equilibrium. 
<5detritus is defined as the difference between detritus in the high (hi) and the low (lo) remineralization cases. Units are mg C m“^ . The STCZ 
is denoted with the crossed line at 40° S and the SAP is denoted with the continuous red line, at about 48° S in the Atlantic Ocean and at 
58° S in the Indian and the Pacific oceans. 


have comparable distributions of surface effective reminer- 
abzation rates; their largest differences amount to about 
30%. GISSEH has lower effective remineralization rates 
than GISSER in the tropics due to its lower SSTs there (about 
10 % less), in the Southern Ocean close to Antarctica, and the 


deep convection regions in the North Atlantic. GISSEH ex- 
hibits higher effective remineralization at the Sub-Antarctic 
Front (SAF, also known as Antarctic Polar Front, approx- 
imately along 50° S, on the northern side of the Antarctic 
Circumpolar Current, ACC) region. At 500 m (Fig. 4, second 
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Fig. 4. Effective remineralization rate at the surface, 500, 1000 and 3000 m for a = 0.1 d ^ (the low remineralization case). The left column 
is for the GISSER model, the middle column for GISSEH and the right column is the difference GISSEH— GISSER. 


row panels) the differences are larger, up to 50 %, for most of 
the ocean. GISSEH mixes down heat more vigorously in the 
deep convection region of the North Atlantic, in the north- 
western Pacific and in the northern Indian Ocean. This ef- 
fect is reversed below 1000 m (Fig. 4e, I) where the GISSER 
ocean is warmer than the GISSEH ocean except in a nar- 
row latitude band adjacent to Antarctica. Therefore effec- 
tive remineralization is 20-60% larger in GISSER than in 
GISSEH. In the abyssal ocean, at 3000 m (Fig. 4g, h), ef- 
fective remineralization is about 30-40% higher in GISSER 
than in GISSEH, except in the narrow latitude band against 
the coast in Antarctica where GISSEH effective remineral- 
ization is about 10% larger than GISSER. Depicted in the 
right-hand column in Fig. 4 are the differences between the 
effective remineralization in GISSEH and GISSER. Away 
from polar latitudes, GISSER is warmer and therefore has 
higher remineralization rates at all depths. Only in the North 
Atlantic subtropical and midlatitudes is GISSEH warmer at 
depth (down to 1000 m at least) and therefore has higher rem- 
ineralization rates. 

The settling speed term in NOBM is parameterized as fol- 
lows: 

w = wsd* (0.451 -F 0.0178 T) (5) 

for ice-free ocean. The rate wsd is 20 md“ * for particulate 
nitrogen/carbon. 


The vertical flux of detritus (normalized to the flux at 
75 m) is very sensitive to the remineralization rate changes. 
Both model flux profiles follow the Martin curve (Fig. 5) 
for the lower remineralization rate cases. Agreement with 
the Martin curve, i.e., the exponential 0.86 form, extends to 
400 m depth for the GISSER model and to 200 m depth for 
GISSEH. Below these depths, however, the models under- 
estimate the Martin curve significantly. This result implies 
that perhaps our models remineralize detritus much faster 
than they allow it to sink. However, as described in Bues- 
seler et al. (2008), there may be large variation in the Martin 
exponent magnitude in the world ocean, at different locations 
and depths. Other studies have also shown that the exponen- 
tial may vary anywhere in between 0.3 and 2 (Kwon et al. 
(2009) and references therein). 

3.1.3 Nutrients 

Similar to detritus, nutrients and in particular nitrates are 
modeled as follows: 

RaD D 

dNt = V{KVN)-V \N - *N[S,7.t; /(NOs); R] -h -h Ad (6) 

with terms analogous as in Eq. (3) and, in addition, C : N 
is the carbon to nitrogen Redfield ratio, is the nitro- 
gen : chlorophyll ratio, pi is the phytoplankton maximum 
growth rate for each species / (i.e., diatoms, chlorophytes, 
cyanobacteria, and coccolithophores) and /(NO3) is the 
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Normalized detritus flux 



Fig. 5. Globally averaged normalized detritus flux (in 
mmolm“^d“*) and comparison to the theoretical Martin 
curve. In blue: GISSER, in red: GISSEH. 


nutrient uptake funetion. The parameter values are the ones 
used in Gregg and Casey (2007); namely, C : N ratio = 
106/16, nitrogen : chlorophyll ratio = 50 (C:N ratio), phy- 
toplankton maximum growth rate for diatoms = 2.00 (in 
units d“^ and at 20 °C), for chlorophytes, diatom growth 
rate = *0.84, for cyanobacteria, diatom growth rate = *0.670 
and for coccolithophores (E. huxleyi only), diatom growth 
rate = *0.755. 

Because detritus at depth is a major source of remineral- 
ized nutrients (Eq. 2), we expect that increased remineraliza- 
tion rate leads to decreases in the amount of detritus at depth, 
and subsequent reductions in the recycled nutrient load of the 
deep waters. Therefore Fig. 6a and b shows lower concentra- 
tions of nitrates at the surface in the deep convection regions 
of the Northern Flemisphere and in the Southern Ocean, ex- 
tending from the coast of Antarctica to the SAF, because of 
mixing with nutrient depleted deep water. Detritus-strapped 
water masses that ventilate the world ocean produce less re- 
cycled nutrients and result in lower nutrient distributions in 
the main upwelling zone in the Southern Ocean. This is a 
robust result in both models. 

Flowever, unlike detritus, north of the STCZ, nutrient con- 
centrations are higher in the high remineralization rate cases. 
This is because these waters originate from the tropical ocean 
upwelled surface and are therefore richer in remineralized 
nutrients. In fact, the increases in remineralized nutrient sup- 
ply at the surface along the STCZ and the South American 
west coast are twice as high as the background values. 

Furthermore, the change in nutrients due to recycling 
changes in the Pacific Ocean is not robust in the two models 
(Fig. 6a, b). GISSER gives more nutrients in the upwelling 
zones at the Equator than GISSEFI. The reason for that is that 


A. Romanou et al.: Ocean biological pump and recycling 

GISSER has much stronger coastal upwelling due to stronger 
winds over the cold tongue region (Romanou et al., 2013), 
which brings more nitrate to the surface. 

As a result, increasing remineralization causes increase of 
nutrients in the equatorial eastern Pacific in GISSER but not 
in GISSEFI where the distributions are already too low. 

Similarly, in the northeastern Pacific, higher remineraliza- 
tion rates in GISSER than in GISSEH, due to higher SSTs in 
this model, result in higher nitrates in the Russell model but 
lower in the HYCOM (Fig. 6a, b). The mixing and the mixed 
layer depth in GISSEH extends deeper than in GISSER and 
the combination of remineralization and local mixing result 
in GISSEH showing a reduction in nutrients as we increase 
the remineralization rate. 

In the subpolar North Atlantic, increasing the remineral- 
ization rate reduces nutrients at the surface because deep 
mixing (through seasonal convection) lowers surface val- 
ues, since less recycled material reaches the deeper ocean 
(Fig. 6a, b). 

The vertical distributions of nitrates show that increasing 
the remineralization rate produces a deeper, less steep nutri- 
cline because more recycled nutrients are generated in the 
upper 1000 m due to more rapid detritus remineralization. 
Vertical mixing of nutrients due to the local nutrient gradi- 
ent at the base of the nutricline is enhanced. The increase at 
the surface occurs mostly in the Southern Ocean (Fig. 6c, d) 
and amounts to about 20 % of the background values. 

Therefore, mesopelagic reductions (at depths of 100- 
1500 m, Fig. 6c, d) are a robust feature in both models, al- 
though in GISSEH the vertical extent is shallower and not 
as pronounced. The largest reductions in nutrients occur in 
the region between the SAF and the STCZ in each model 
where the nutricline slopes upward. The GISSEH subsurface 
nutrient reductions are more pronounced, because in GIS- 
SEH upwelling extends deeper than in GISSER. The region 
of the largest sensitivity to changes in the remineralization 
rate is the mesopelagic “twilighf’ zone as has been previ- 
ously pointed out by Buesseler et al. (2008). 

These results are robust across the models particularly in 
the 30°-50° S belt. Thus, in this area nutrient distributions 
are primarily controlled by sensitivities to biological param- 
eterizations rather than differences between the two physical 
models, which mostly are due to physical parameterizations. 

3.1.4 Total chlorophyll, primary production and deep 
carbon export 

Chlorophyll distributions are determined by light and nutri- 
ent availability and the temperature of seawater. GISSEH 
predicts higher chlorophyll concentrations than GISSER 
(Fig. 7) in two main regions: (1) the tropical upwelling re- 
gions, where water originates from deeper, generally richer 
in nitrates, layers and (2) the Southern Ocean (south of 
50° S) where iron concentrations are higher. However, GIS- 
SEH has lower chlorophyll concentrations in the subtropical 
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Fig. 6. (a), (b) Surface distributions of ^Nitrate in (a) the GISSER and (b) the GISSEH models at equilibrium. SNitrate is defined as the 
difference between nitrate in the high and the low remineralization cases. Units are in mmol N m“^ . (c), (d) Zonally averaged vertical sections 
of iSnitrate defined as the difference in nitrate between high and low remineralization rates in the (a) GISSER and (b) GISSEH models. Units 
are in mmolN m“^. 


gyres because of steeper nutriclines (as well as thermo- 
clines/haloclines), which prohibit vertical mixing with the 
deeper layers, resulting in unfavorable conditions for phy- 
toplankton growth. 

In the Southern Ocean, in the region between the SAP 
and the STCZ, increasing the remineralization rate leads to 
increased nutrient supply, which in turn results in 30-50 % 
increase in total chlorophyll in both ocean models (Fig. 7). 
South of the SAP, however, increased recycling of nutrients 
reduces nutrient concentrations (Fig. 6c, d) which in turn 
leads to lower surface chlorophyll concentrations. For similar 
reasons, in the high latitude Atlantic ocean, phytoplankton is 
also reduced (see Figs. 6a, b and 7e, f). 

In the tropicaPsubtropical ocean (40° S-50° N), chloro- 
phyll changes with remineralization are robust: both mod- 
els agree quantitatively in the magnitude of chlorophyll as 
well as the fact that increases in remineralization rate lead 
to increased chlorophyll (Fig. 8a). In the subpolar regions, 
however, the two physical models, with different high lati- 
tude treatment of deep convection, ice formation and brine 
rejection, as we have seen in the Model Description section, 
lead to quite different distributions of chlorophyll. It is found 
that at high latitudes, the chlorophyll sensitivity to the physi- 


cal parameterizations is larger than the sensitivity to changes 
in the remineralization rate. Overall, increasing the reminer- 
alization rate tends to increase surface chlorophyll, except in 
the Southern Ocean south of the SAF, in the region 40-70° S, 
because of lower nutrient supply there (Fig. 6a, b) from the 
upwelling deep water. 

Primary production (PP) is very different in the two mod- 
els. HYCOM has consistently lower PP than the Russell 
model, because of the much shallower mixed layer depths. 
However, the response to increased recycling rates is robust 
in both models, with primary production increasing with the 
recycling rate in all the regions (Fig. 8b). These increases are 
important in the Northern Hemisphere’s frontal regions 30- 
60° N, in the subtropical gyres, because of the more readily 
available remineralized nutrients there. Little change is found 
in the Southern Hemisphere’s high latitudes, south of 40° S, 
where changing the remineralization rate does not impact PP 
(Fig. 8b). Diatoms are less influenced by the remineraliza- 
tion rate changes in both models. Cyanobacteria and coc- 
colithophores are the most sensitive particularly in the tropi- 
cal/subtropical regions and chlorophytes mostly in the North- 
ern Hemisphere’s subtropical front region. The reason for 
the disparate response of the chlorophyll species to recycling 
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(a) R-lo: Tchl at the surface (b) H-lo; Tchl at the surface 
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Fig. 7. Surface distributions of chlorophyll in (a) GISSER and (b) GISSEH with low remineralization rate, (c) GISSER and (d) GISSEH 
with high remineralization rate. Units are in mg C m“^ . (e), (f) Differences in surface chlorophyll between same model runs but with varying 
remineralization rate. Units are in mg C m“^ . 


changes is that growth, which is the process that is affected 
most by remineralization changes, depends on the geograph- 
ical distributions of species and nutrient. Diatoms, which are 
mostly found in the high latitude oceans in our models, ben- 
efit less from the increase in nutrients. Coccolithophores and 
cyanobacteria, which are more abundant than diatoms par- 
ticularly in low and midlatitudes, benefit more from the in- 
creases in nutrients when the remineralization rate increases. 

Carbon export, i.e., the fiux of particulate organic material 
across 75 m depth, illustrates the strength of the biological 
pump. The carbon export is very sensitive to both physical 
and biological parameterizations (Fig. 9a, b). Altering the re- 
cycling rate in each model leads to changes in carbon export. 
These changes are more pronounced in the GISSER model, 
however, the relative changes between the high and the low 


remineralization rate runs are of the same order of magnitude 
in both models (Fig. 9a). We find that an 80 % increase in the 
remineralization rate leads to about a 50-70 % reduction in 
the carbon export. Carbon export decreases as the remineral- 
ization rate increases, which is due to the reduction of partic- 
ulate organic matter concentration at depth as more of it rem- 
ineralizes near the surface. The largest differences between 
the high and low remineralization runs occur at the frontal 
regions of the subtropical gyres and at the Equator, which are 
the same regions where detritus differences are found to be 
largest (Fig. 3). In these regions, the two physical model re- 
sponses are closer than the differences in response due to the 
biogeochemistry. This becomes apparent from the compari- 
son of the vertical fluxes across the 75 m isobath (Fig. 9b). 
These vertical fluxes are associated with processes such as 
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(a) R-H: zonally avgd tctil at surf 


(b) R-H: zonally averaged PP 


) PP VS remineraiiz. rate 


(b) cexp vs remineraiiz. rate 



Fig. 8. Zonally averaged (a) total chlorophyll at the surface 
(mgchlm”^) and (b) primary productivity vertically integrated 
(mgCm“^d“*). Continuous lines are for high remineralization 
rate and dashed for low remineralization rate. Blue lines are for 
GISSER and red lines for GISSEH. 


Vertical fluxes across the 75m Isobath 



Fig. 9. (a) Zonally integrated carbon export at 75 m in Pg Cyr“ * . 
Color scheme is same as in previous figures, (b) Zonal mean vertical 
fluxes across 75 m depth for the Russell model: carbon export (R-lo, 
continuous line), vertical advection flux (multiplied by 10, dashed 
line), vertical flux associated with the Gent-McWilliams scheme 
(multiplied by 10, solid line with crosses - not well visible on this 
diagram as it coincides with the axis) and vertical flux associated 
with vertical mixing (multiplied by 10, solid-circle line). All fluxes 
are positive down. 


vertical mixing (KPP), vertical advection (which is also sen- 
sitive to the advection scheme) and the Gent-McWilliams 
mesoscale eddy mixing. It is found that in the Russell ocean 
model, all physical circulation related vertical fluxes are 1-2 
orders of magnitude smaller than the carbon export at that 
level. 

Additional model runs with different remineralization 
rates were performed to examine the sensitivity of the ocean 
biological pump to a wider range of values for this param- 
eter. The results show that both models (GISSER and GIS- 
SEH) exhibit the same sign of the sensitivity, although the 
magnitude is different (see Fig. 10a, b) for the primary pro- 
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Fig. 10. (a) Vertically averaged PP as a function of remineralization 
rate. The vertical axis is global integral PP in Pg Cyr“ ^ and the hor- 
izontal axis is remineralization rate per day. Blue line is for GISSER 
and the red line for GISSEH. (b) Globally integrated carbon export 
across 75 m (cexp, in PgCyr“^) as a function of the remineraliza- 
tion rate (d“ * ), (c) globally integrated efficiency of the biological 
pump (ratio of the carbon export at 75 m and primary productivity). 
The bars are estimates from the models and the symbols represent 
analytic solutions obtained using Eq. (7) for each model, PP from 
diatoms, PP from chlorophytes, PP from cyanobacteria and PP from 
coccolithophores. 


duction and the carbon export: increasing the remineraliza- 
tion rate leads to increased primary production and decreased 
carbon export. Primary production responses (Fig. 10a) show 
a tendency towards saturation at remineralization rate val- 
ues of about 0.3 and higher, whereas carbon export does 
not (Fig. 10b). The reason for this is that there is a limit 
of nutrient remineralization increase, beyond which, limita- 
tions of other nutrients (iron and silica) start playing a role. 
Carbon export however continues to decrease as remineral- 
ization rate is increased, because particulate organic carbon 
continues to decrease. 

Estimates of the f ratio (or carbon export efflciency, i.e., ra- 
tio of export at 75 m to the vertically averaged PP) are shown 
in Fig. 10c: the carbon export efficiency in both models de- 
creases as we increase the remineralization rate. Model car- 
bon export efficiency estimates (0. 1-0.7) fall within observa- 
tional estimates as reported in Laws et al. (2000) and Henson 
etal. (2012). 

Both models show similar rates of decrease of the f ratio 
with increasing remineralization although the magnitude of 
the f ratio itself is very different in the two models. To fur- 
ther investigate the sensitivity of the carbon export efficiency 
to model parameterizations, both physical and biological, we 
used parametric curve fitting that describes the functional 
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relationship between the f ratio and the remineralization in 
the two models. Initially, two parametric equations were ob- 
tained, one for each physical model. The two sets of coeffi- 
cients were then expressed in terms of a single metric of the 
physical models, here chosen to be the model mixed layer 
depth normalized by deBoyer Montegut et al. (2004) ob- 
served mixed layer depth. The choice of mixed layer depth in 
expressing the dependence of export efficiency on remineral- 
ization is probably not unique, but was dictated by the fact 
that MLD is an expression of the physical model’s internal 
dynamics and ability to redistribute tracers. Different forms 
of the mixed layer depth such as the globally averaged mixed 
layer depth, the North Atlantic mixed layer depth, which is 
associated with the strength of deep convection in that region, 
and the Southern Ocean mixed layer depth, which is the re- 
gion were remineralized nutrients resurface, were tried. The 
metric that provided us with the best fit in the least square 
sense to both model estimates of the f ratio at all the rem- 
ineralization rate values was the normalized Southern Ocean 
averaged mixed layer depth (see Eq. 7). 

Equation (7) describes how the biological pump efficiency, 
Abp, changes with the remineralization rate a: 

= 0.017 rMLDso - 0.017 rMLDso « -I- 0.06 rMLDso (7) 

where rMLD is the Southern Ocean normalized mixed layer 
depth (annual mean, regional average 70°-40° S normalized 
by observations) in each model. Equation (7) implies that al- 
though there is a tight relationship between the remineraliza- 
tion rate and the carbon export efficiency, this is modulated 
by how well the physical model depicts mixing in the South- 
ern Ocean. It would be interesting to see if this or similar 
functional dependencies of the biological pump efficiency to 
physical and biological parameterizations are valid in other 
biogeochemical models. 


3.2 Sensitivity of the gas exchange pump 
3.2.1 Die and pC02 

Changes in remineralization impact CO 2 gas exchange be- 
tween the ocean and the atmosphere through changes to 
surface DIG and / 7 CO 2 . Similar to nitrates (Fig. 6a, b), in 
both models increasing the remineralization rate leads to in- 
creases of DIG in the region between the SAF and the STCZ 
(Fig. 1 1), again due to the fact that carboclines get steeper (as 
do nutriclines), and therefore mixing with underlying waters 
is suppressed. In the Southern Ocean divergence zone, DIG 
increases are small as a percentage of the total DIG (about 
10%, Fig. 11) but they are of the same order of magnitude 
as the local meridional gradient. South of the SAF, towards 
Antarctica, there are lower concentrations of DIG in the high 
remineralization case, due to the upwelling of deep carbon- 
deficient water. As in nitrates, the increase of surface DIG 
with the remineralization rate in the tropics is found only 
in GISSER, and is due to the more nutrient rich waters in 
GISSER than in GISSEH. 

In the subtropicaFequatorial Pacific the two models dis- 
agree. GISSER predicts increased DIG at the surface in the 
runs where recycling is increased while GISSEH shows a 
weak reduction. The disparity is attributed to two competing 
processes that set the DIG distribution at the surface when we 
increase the remineralization rate: remineralization sources 
vs. vertical mixing and sinking. As the surface remineralized 
DIG distribution increases as we increase the rate of reminer- 
alization, the subsurface profile of DIG changes as well due 
to sinking. As we increase the remineralization rate, subsur- 
face DIG decreases. Vertical mixing (or upwelling) will then 
dilute the surface distributions with deeper and lower distri- 
butions of DIG. Here, GISSER has higher temperatures than 
GISSEH and therefore more remineralized DIG at the sur- 
face. As we increase the remineralization rate, initially both 
models will have more remineralized DIG at the surface; but 
the GISSER depth profile of remineralized DIG will be more 
surface enhanced and GISSEH will have more remineralized 
GISSER below the surface. Vertical mixing comes into play 
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(a) R: zonally avgd diff DIG section (Rhi-Rlo) (b) H: zonally avgd diff DIG section (Hhi-HIo) 



Fig. 12. Differences in the DIC zonally averaged vertical sections between same model mns with varying remineralization rates. 


then, as GISSEH has deeper mixed layers (see Fig. 7 in Ro- 
manou et al. (2013)) and therefore mixing with subsurface 
DIC is such that GISSER ends up with more DIC whereas 
GISSEH ends up with less DIC. 

Similar to nitrates (Fig. 6c, d), increasing the reminer- 
alization rate leads to decreased concentrations of DIC at 
mesopelagic depths (100-2000 m. Fig. 12). As detritus is 
remineralized back to inorganic carbon, less of it descends 
through the water column, reducing the source of reminer- 
alized carbon at depth. This is a robust behavior, evident in 
both physical models. In GISSER, increasing the remineral- 
ization rates deepens the carbocline and makes it more dif- 
fuse. This effect extends down to 1500 m depth. In GISSEH, 
there is a similar effect but it only extends down to 500 m. 

In both models, as we increase remineralization rates, 
the carbocline gradient at mesopelagic depths decreases and 
therefore more DIC can reach the surface through vertical 
mixing, which explains the increased DIC distributions there. 
This is mostly evident north of 50° S (i.e., north of the SAF). 
The changes in partial pressure of CO 2 at the surface of the 
ocean (pC02) mirror the changes in the surface DIC with the 
remineralization rate. 

3.2.2 Equilibrium flux of CO 2 

In the model runs described here, changes in the surface 
/ 2 CO 2 depend on changes in DIC due to varying the rem- 
ineralization rate. Surface alkalinity is modeled as a constant 
function of surface salinity and is also affecting the flux. 
However, since we do not simulate the carbonate pump in the 
model we do not account for changes in alkalinity at depth. 
Since both models are at equilibrium, T, S and surface winds 
do not change within runs of the same model as we vary the 
remineralization rate. Therefore we are able to examine how 


much gas exchange depends on changes in DIC and the bi- 
ological pump. Note that in the following, the flux sign con- 
vention is positive down, i.e., positive air-sea flux denotes 
oceanic uptake of CO 2 . 

The Southern Ocean sink (Fig. 13) is substantially reduced 
in the fast recycling case. This reduction occurs mainly be- 
cause the uptake region (between 40° and 60° S) is signif- 
icantly reduced while although, the source region south of 
40° S is also reduced, this reduction in outgassing is not suffi- 
cient to outweigh the reduction of the actual sink more to the 
north. Carbon uptake is reduced because increased reminer- 
alization leads to higher DIC and PCO 2 at the surface at the 
divergence zone of the Southern Ocean. South of the SAF, 
outgassing is reduced because although surface DIC (and 
PCO 2 ) is reduced, due to carbon deflcient NADW, the im- 
balance with the atmosphere still favors outgassing, but less 
so than in the low remineralization case. At the same time, 
tropical outgassing is similar in the two models, although 
slightly enhanced in GISSER due to more upwelling DIC. 
The North Atlantic sink strengthens because surface DIC de- 
creases as increased remineralization leads to less carbon at 
depth. Therefore more CO 2 is sequestered into the ocean. On 
the other hand, the North Paciflc sink of CO 2 gets stronger 
in the GISSER model when we increase the remineralization 
rate, but the opposite happens in the GISSEH model. 

As a global average, the outgassing flux increases as we 
increase the remineralization (region Gib in Fig. 14). This 
increased outgassing occurs mostly in the Southern Ocean 
(ANT in Fig. 14) and in the equatorial regions (EIN, EPA 
and EAT in Fig. 14). The North Atlantic sink increases in the 
GISSER model but slightly decreases in GISSEH because 
mixing in GISSEH brings much less DIC near the surface. 
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Fig. 13. Equilibrium flux (mol CO2 m“^ yr~'). Red denotes uptake, blue denotes outgassing. In unshaded regions the flux is low. These are 
results from the two models GISSER and GISSEH using different remineralization scales. Hi refers to high remineralization rate (i.e., fast 
recycling of nutrients) and lo refers to low remineralization rate. 


4 Discussion 

In this paper we assess the sensitivity of the oceanie carbon 
cycle and the air-sea CO 2 exchange to changes in the recy- 
cling rates within the ocean. We have also compared these 
changes with model uncertainties due to varying parameter- 
izations of physical mechanisms such as vertical advection 
due to large-scale flow convergence/divergence, mesoscale 
eddy stirring and turbulent mixing. For this we ran the GISS 
climate model coupled to two different ocean models and 
we varied the remineralization rate gradually from 0.01 to 
0.5 d“*. The two ocean models used in this study differ in 
the vertical coordinate system and therefore simulate diapy- 
cnal and isopycnal processes differently. 

The main result obtained is that carbon cycle sensitivities 
to the biological pump parameterizations are of the same or- 
der of magnitude or larger than the sensitivities due to phys- 
ical model parameterizations. Consequently, equal attention 
should be placed in assessing the fidelity of parameteriza- 
tions such as remineralization scales, sinking rates and oth- 
ers that may affect the biological pump as to those that aim 
to improve the physical ocean skill. Specifically, organic par- 
ticulate sinking reductions due to increased recycling in the 
epi- and mesopelagic ocean are larger than uncertainties in 
parameterizations of physical processes, particularly in sub- 
tropical frontal regions in the Northern and Southern Hemi- 


spheres. This leads to decreased nutrient remineralization at 
these depths and nutrient depletion of deep water masses 
such as the NADW. Primary production increases as we in- 
crease the remineralization rate while carbon export at 75 m 
decreases. CO 2 outgassing is enhanced as we increase the 
remineralization rate although CO 2 uptake increases in the 
SAF region. 

In this study, the Southern Ocean CO 2 sink emerges as par- 
ticularly sensitive to the parameterizations of the biological 
pump, and this is a robust feature in both models. In par- 
ticular, waters that upwell in the Southern Ocean and which 
originally mixed and subducted in the Northern Hemisphere 
(NADW) carry less nutrients and DIG in the runs where re- 
generation rates are higher. The Southern Ocean between 
APF and SAF supports less chlorophyll and outgasses less 
CO 2 to the atmosphere in these runs as well. North of the 
SAF in the Southern Ocean, between SAF and the subtrop- 
ical front (STF), surface DIG is increased and the sink re- 
gion there is reduced. Ghanges in the biological pump in the 
Southern Ocean’s and the Northern Hemisphere’s frontal re- 
gions lead to substantial changes in the surface DIG distri- 
butions and therefore in the air-sea exchanges. This result is 
robust in both ocean models and therefore independent of a 
model’s physical parameterizations. 

The sub-Antarctic zone (SAZ, 35-50° S; region between 
SAF and STF), is an area of mode water formation and the 
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Table 1. Correlation coefficients of the spatial patterns between model experiments and observations. Values in parentheses are below the 
5 % confidence level. 


Model exp 

tchl 

diat 

chlo 

cyan 

cocc 

nitr 

sili 

iron 

die 

GISSER-high rem 

0.64 

0.78 

0.86 

0.77 

(0.39) 

0.98 

(0.36) 

(0.25) 

0.84 

GISSER-low rem 

0.65 

0.74 

0.83 

0.83 

0.68 

0.99 

(0.45) 

(0.30) 

0.83 

GISSEH-high rem 

(0.52) 

0.75 

0.86 

0.59 

0.74 

0.97 

0.61 

(0.18) 

(0.50) 

GISSEH-low rem 

(0.42) 

0.71 

0.85 

0.60 

0.74 

0.98 

0.61 

(0.18) 

(0.50) 



Fig. 14. Regional fluxes (in mmol C m~^ yr~'). The values are 
scaled to take into account the different areas of the regions. ANT: 
around Antarctica, SIN: southern Indian Ocean, SPA: southern Pa- 
cific Ocean, SAT: southern Atlantic Ocean, EIN: equatorial Indian 
Ocean, EPA: eastern Pacific, EAT: eastern Atlantic, NIN: north- 
ern Indian, NCP: north-central Pacific, NCA: north-central Atlantic, 
NPA: North Pacific, NAT: North Atlantic and Gib: global ocean. 


most vigorous atmospheric sink of CO 2 (Metzl et al., 1999; 
Sabine et al., 2004; Fletcher et al., 2006; Takahashi et al., 
2009). The SAZ is also the location of enhanced phytoplank- 
ton blooms as opposed to the HNLC regions in the rest of the 
Southern Ocean (Lourantou and Metzl, 2011). The South- 
ern Ocean, southward of 30° S, is a major uptake region for 
atmospheric CO 2 , in atmospheric and ocean inversion mod- 
els (Friedlingstein and co authors, 2006; Gruber et al., 2009) 
as well as in observational studies (Metzl et al., 1999; Taka- 
hashi and CO authors, 2002; Takahashi et al., 2009) although 
the magnitude of this uptake is argued (Gruber et al., 2009; 
Denman and Brasseur, 2007) . 

Mixed layer depths in this region where mode waters 
form, e.g., the Sub-Antarctic Mode Water (SAMW) and the 
Antarctic Intermediate Water (AAIW), increase particularly 
along the SAF, but the uncertainty in the observations is large 
since this region is severely undersampled, especially dur- 
ing austral winter. Our models show that increasing reminer- 
alization rate leads to modest decreases in pC 02 of about 


6-10%, and the Southern Ocean CO 2 flux sensitivity to the 
biological pump parameterizations is greater than model un- 
certainty due to physical processes. 

As in previous studies (Kriest and Oschlies, 2011; Mari- 
nov et al., 2008), we conclude that in the Southern Ocean, 
the efficiency of the biological pump determines the atmo- 
spheric pC02 distributions. Surface nutrients are not be- 
ing efficiently utilized in the Southern Ocean and this leads 
to deep water masses (which have been ventilated in the 
Southern Ocean) that have more preformed nutrients (i.e., 
nutrients that have not been biologically utilized, originate 
from the surface and are redistributed through mixing, “un- 
remineralized nutrients”) than NADW. As the deep water 
masses, rich in preformed nutrients and ventilated in the 
Southern Ocean, upwell at low latitudes, they promote in- 
creased bio-utilization of surface nutrients and more carbon 
fixation at the surface. Flowever, increasing remineralization 
depth (or decreasing remineralization rate) depletes surface 
nutrients in the NADW more than in the Southern Ocean, be- 
cause in the latter the nutricline is less steep. Consequently, 
deep water masses ventilated from the Southern Ocean are 
not as efficient as NADW at storing atmospheric CO 2 when 
the remineralization rate decreases, a result in agreement 
with Kwon et al. (2009). 

Carbon export efficiency is determined by the fraction 
of carbon fixed by primary production that is transported 
through recycling and sinking to the deep ocean. This carbon 
flux is important because it removes carbon that is in contact 
with the atmosphere and sequesters it to greater depths and 
therefore increases the ocean’s ability to uptake more CO 2 . 
The global mean export ratio is estimated to be between 10 
and 40% (Flenson et al., 2012) and both GISS models lie 
within this range. 

The globally integrated carbon export is found to be about 
5±2.2PgCyr“^ (Flenson et al., 2012), 6PgCyr“* (Martin 
et al., 1987) or even 10PgCyr“* (Dunne et al., 2007). We 
And carbon export ranging from 1 to 7 Pg C yr“ * according to 
remineralization rates and model physics. GISSEFl has typ- 
ically less carbon export across the 75 m depth and as we 
increase remineralization rate, the carbon export decreases. 
The carbon export efficiency (i.e., the ratio of carbon export 
to primary productivity, Sarmiento et al. (2004)) is about 0.5 
for the low remineralization rate in GISSER and drops to 
0.17 for the higher rates we simulated. GISSEFl is generally 
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half as efficient in exporting carbon from the euphotic zone. 
Our remineralization profiles agree with Martin et al. (1987), 
who found that about half of the carbon exported from the eu- 
photic zone towards greater depths is regenerated at depths 
shallower than 300 m; whereas two-thirds is regenerated in 
the upper 500 m and 90% in the upper 1500 m. This very 
fast recycling in the euphotic zone implies that the majority 
of the carbon exported from the surface will be again in con- 
tact with the atmosphere within tens of years or less (Martin 
etal., 1987). 

In addition to carbon, nutrient distributions are also sen- 
sitive to changes in remineralization rate, particularly in the 
SAZ (Pollard et al., 2002; Lovenduski and Gruber, 2005), 
the formation region of the SAMW, which supplies most 
of the world ocean surface water with nutrients (Sarmiento 
et al., 2003). This area has been shown to be very sensitive to 
changes in the Atlantic Meridional Overturning Circulation 
(AMOC; Schmittner (2005)) and in our results is shown to 
be very sensitive to biogeochemical changes as well. The fact 
that HYCOM has shallower penetration of ventilated NADW 
than the Russell model affects remineralization levels, acting 
to reduce the remineralization depth. 

Finally, although the scope of this paper is not to evalu- 
ate the quality of the model predictions against observations 
(since the model reflects preindustrial equilibrium) a brief 
comparison to observed fields is given in Table 1. In particu- 
lar, regional correlations of the annual mean surface distribu- 
tions in total chlorophyll, diatoms, chlorophytes, cyanobac- 
teria, coccolithophores as well as nutrient fields (nitrate, sil- 
icate and iron) and DlC are shown in Table 1. GISSER 
yields equal or better regional correlations with observa- 
tions in the low remineralization case for all these fields 
whereas GISSEFI distributions seem to be much less af- 
fected by the change in the remineralization rate. While this 
is not a rigorous model evaluation, it is evident that the bi- 
ological sensitivity of the two models is different with the 
GISSER model being more transparent to these changes. 
The observational data for the estimation of correlation co- 
efficients come from the data set described in Gregg and 
Casey (2007), Sect. 2.3. This data set includes 469 surface 
distributions of phytoplankton abundances and is available at 
http://polar.gsfc.nasa.gov/research/oceanbiology/index.php. 


5 Conclusions 

Simulations of the natural carbon cycle using ocean mod- 
els that parameterize differently physical processes such as 
convection, diffusion and advection of thermodynamic prop- 
erties and tracers, can be very instmctive. They help set the 
range of uncertainty in the oceanic CO 2 uptake, nutrient and 
carbon distributions, both in global but also in regional terms. 
Varying a biogeochemical parameterization of the ocean bi- 
ological pump, such as the recycling rate, can help us assess 
how important biogeochemical processes are in controlling 
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the air-sea flux of CO 2 , the regional distributions of CO 2 
sources and sinks as well as the carbon export efficiency, 
compared to physical model uncertainties. 

We found that uncertainties due to biology become as im- 
portant and occasionally more important than uncertainties 
in physics and this may suggest that improvements to phys- 
ical models and biogeochemical models should be carried 
out in parallel. In our study, the Southern Ocean emerges as 
particularly sensitive to changes in the recycling rates and 
less so to physical model differences. The magnitude of the 
global biological pump depends on the effectiveness of mix- 
ing in the Southern Ocean but its sensitivity to recycling is re- 
markably robust, i.e., independent of physical parameteriza- 
tion choice. Similar analyses to more ocean/biogeochemical 
models could further generalize these results. 
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